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ABSTRACT 

To overcome the limitations of classical propeller theory, a computer 
program, NASPROP-E, has been developed which solves for the flow field sur- 
rounding a multi-bladed propeller and axisymmetric nacelle combination using a 
finite-difference method. The governing equations are the three-dimensional 
unsteady Euler equations written in a cylindrical coordinate system. They are 
marched in time until a steady state solution is obtained. The Euler equa- 
tions require no special treatment to model the blade work vorticity. The 
equations are solved using an implicit approximate factorization method. 
Numerical results are presented which have greatly increased the understanding 
of high speed propeller flow fields. Numerical results for swirl angle down- 
stream of the propeller and propeller power coefficient are higher than exper- 
imental results. The radial variation of swirl angle, however, is in reason- 
able agreement with the experimental results. The predicted variation of 
power coefficient with blade angle agrees very well with data. 
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** Research Scientist; member AIAA. 

*** Chief, Applied Computational Aerodynamics Branch; Associate Fellow, AIAA. 
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NOMENCLATURE 


A 

b 

B 

c 

C 

CN 

D 

e 

E 

F 

G 

H 

H r 

( i z > V»\ ) 

J 

k-j , i =0 , 3 
LEA 
n 
P 


Jacobian matrix, 3E/3Q 
blade width or chord 
Jacobian matrix, 3 F/3 Q 
speed of sound, (-rp/p)^'^ 

Jacobian matrix, 3G/3Q 

Courant number (see Eq. (7)) 

diameter of propeller 

total energy per unit volume 

vector of flux quantities in ^-direction 

vector of flux quantities inn-direction 

vector of flux quantities in c-direction 

vector of source terms created by cylindrical 

generalized coordinate transformation 

rothalpy 

unit normal vector of cylindrical coordinate system 

Jacobian or advance ratio, (U^/nD) 

generalized matrices 

leading edge alignment 

rotational speed, revolutions per second 

static pressure 

2 

dynamic pressure, O.Bp^v^ 


Qi, 1=1,5 

Q 

R 

(t,z,r,0 ) 

u 

U 

v 

V 


components of Q 

vector dependent variable of integration 
blade tip radius 

independent variables, cylindrical coordinates 
physical velocity in z-direction 
contravariant velocity defined in Eq. (2) 
physical velocity in or r-direction 
contravariant velocity defined in Eq. (2) 


V velocity vector, ui z + vi r + wi^ 

7 vector form of the contravariant velocity 

w physical velocity in or ^-direction 

W contravariant velocity defined in Eq. (2) 


33/4 

blade angle at r/R = 0.75 


Y 

ratio of specific heats 


AS 

blade twist 


AP 

incremental pressure 


AS, An, AC 

computational mesh spacing 


P 

density 


0 i, 1=1,5 

eigenvalues of gas-dynamic equations (see Eq. (6)) 

(t ,5 ,n ,C ) 

transformed independent variables 

(Eq. (1)) 

e e 

explicit smoothing coefficient in 

implicit algorithm 

e i 

implicit smoothing coefficient in 

implicit algorithm 

W 

angular velocity 



Subscripts 

i integer mesh point location in ^-direction 
j integer mesh point location in n-direction 
k integer mesh point location in c-direction 
00 free stream conditions 
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INTRODUCTION 


The continuing quest for improved propulsive efficiency for subsonic air- 
craft has caused a growing interest in the propeller as an alternative to the 
turbofan engine at cruise Mach numbers up to 0.8. Numerous mission studies 
have shown that advanced turboprop powered aircraft offer a potential 15 to 25 
percent trip fuel savings over comparable technology high bypass ratio turbo- 
fan powered aircraft at Mach 0.8. Aerodynamic design of propellers for high 
subsonic flight speeds results in a large number of highly swept blades and a 
nacelle and spinner designed to have favorable interference effects on the 
propeller. 

Interest from a computational point of view was stimulated by the need for 
understanding the transonic flow phenomena around propellers and the potential 
for a more efficient computer generated design. Over the years, computational 
procedures for solving fluid flow problems provided an inexpensive but accu- 
rate means of determining the aerodynamic characteristics of complex configu- 
rations. In addition, they have provided the designer with an effective tool 
for maximizing aerodynamic efficiency without the expense of actually building 
and testing numerous designs. Finally, computational methods have often been 
capable of providing information not readily obtainable from experiments. 

The advanced high speed propeller model is a good example of the type of 
configuration for which it is difficult to experimentally obtain aerodynamic 
information needed for performance analyses and design. The blades are virtu- 
ally impossible to adequately instrument because of the high rotational veloc- 
ities, and their relatively small thickness. Also, details of the surrounding 
flow field can only be obtained by performing costly and time consuming flow 
field surveys. Computationally, however, the entire flow field can be deter- 
mined from a single solution of the governing equations including near and 
far-field effects which can be used in acoustic analysis programs and blade 
surface pressure distributions which can be used in structural and aerodynamic 
design analysis programs. 

The theoretical development and subsequent numerical solutions described 
herein are concerned with simulating the inviscid flow about a wind tunnel 
model which consisted of an eight-bladed propeller and spinner with an axisym- 
metric nacelle (instead of the conventional three-dimensional nacelle with 
inlet). Such a configuration thus requires that only the flow about a single 
blade be computed because of periodicity. 


DERIVATION OF GOVERNING EQUATIONS 

To enhance numerical accuracy and efficiency, a nonorthogonal coordinate 
transformation of the governing equations in a particular base coordinate 
system is employed. This maps the surface of the nacelle and both sides of 
the blade onto constant coordinate surfaces which facilitates the application 
of boundary conditions and permits grid point clustering at the body where 
gradients of the dependent variables are expected to undergo rapid changes. 
Use of such transformations, in addition, permits utilization of uniform dis- 
cretization formulas and well-ordered interior grid point solution algo- 
rithms. Under this transformation the equations can still be written in 
conservation-law form to take advantage of the shock capturing properties. 
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The basic orthogonal coordinate system utilized is cylindrical with z 
oriented along the rotational axis, r extending radially outward from the 
z-axis, and $ the meridional angle measured from a vertical plane 
(see fig. 1). It should be reiterated that for this study, only the flow 
between two of the blades is computed, i.e., between the pressure side of one 
blade and the suction side of the next blade, because of periodicity. The 
cylindrical coordinates are transformed to align the blade and nacelle sur- 
faces with various computational planes according to the following: 

t = t 

l = C(t,z,r, 4 >) (1) 

n = n(t,z,r,<t>) 
c = t (t, z,r,4>) 

* 

This generalized nonorthogonal coordinate transformation maps the spinner 
and nacelle onto a constant n-plane and each side of the blade, i.e., the 
suction and pressure sides, onto parts of a constant c-plane. The remaining 
parts of the constant c -planes are periodic surfaces. The radial far-stream 
and outflow boundaries are situated far enough from the prop-fan to minimize 
the reflection of waves. 

The governing partial differential equations in weak conservation-law form 
for cylindrical coordinates under the assumptions of inviscid compressible 
flow and a perfect, non-heat conducting gas for the transformation given in 
Eq. (1) are: 


where 



Q + E + F + G + H = 0 
x £ n c 
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and 
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U, V, and W are the contravariant velocities written without metric normaliza- 
tion. J is the transformation Jacobian and is defined below. Use of the 
Euler equations in conservation laws form guarantees the accurate calculation 
of the shock waves occurring in the flow field. 

In the conservative variables of Eq. (2(a)), the pressure p is nondimen 
sionalized by p^, the density p by p^, and the cylindrical velocity 
components u, v, and w by a to //7 where a^ is the free-stream speed of 

p 

sound (a^ = yp^/p^) and y is the ratio of specific heats. Other 
quantities made dimensionless are the time t by a 00 /(D> / 7) and the angular 
velocity u by D/T/a^. The pressure, density, and velocity components are 

related to the total energy per unit volume e by the following equation for an 
ideal gas: 


e = p/ (r — 1 ) + p(u 2 + v 2 + w 2 )/2 (3) 

The metrics of Eq. (2(b)) are obtained by the chain rule expansion of 
z^, r^, etc. and solved for S z , etc. to yield the following expressions: 
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n 

z 



c z <Vn 
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(r^z c 
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where 


I = 1/J = 




r r 4> z 
K n c 


) r z r 
C n c 


- z. 


) r 
n c 


Vn< 


- <t> _r z 

5n ? 


The quantity <t> T is the rotational or angular velocity of the propeller. 

In general, the metrics of Eq. (4) are not known analytically and must be 
determined numerically at each step of the integration procedure. To accom- 
plish this, second order central-difference formulas are used at interior 
points and three-point one-sided formulas are used at the boundaries. 

The implicit algorithm to be discussed below requires the analytic deter- 
mination of the Jacobians A = 8E/3Q, B = 3F/3Q, and C = 3G/3Q. These Jacobian 
matrices result from the time linearization of E, F, and G in deriving the 
numerical algorithm. The Jacobian matrices can be written in qeneral as 
follows: 


K o 

1 

1 1(1 

K 2 


1 K 3 1.1 

! ' 


1 Ko +7 + K 1 (2 -*>^ 



r 

- K 1 (t ' l > \ + K 2 

q 2 

q l 

1 1 

1 q 4 K 3 q 2 1 


K lqp- K ^-»qf 

Kq + V ♦ K 2 (2 - y) 


q 3 

q l 

q 4 ^3 q 3 

IV- 1 ) 

r V- V - q j V 

1 q 4 ^3 q ? 

| h^-T-li-VTy 

q 4 K 3 

K rq l - 7 1 (T - 11 

q 3 

q l 

r " "1 

1 — ^3 q 4 1 ^3 

K 0 *v +T i(2-Y)^ jr 1 *— A > 





r „ 1 

. * i V + ( T . 1 ) v 2 

1 -(t - 1) - K 1 

q i - 

-(t - 1 ) V ♦ K. 

q 5 

T TT~ 

k q + yV 

T q 1 

1 q l q l 

2 

l 

q l 

1 q 1 r q 1 0 

1 


[ 

- K 2 V 2 


LV(T.-l).* | 

1 r ‘ 1 J 


where V = (k^ + k 2 q 3 + -pr^ q 4 ) /Qi and V 2 = (q| + + qjj ) /qf * 

The terms q-j, i = 1,5 are the components of the dependent variable of inte- 
gration Q in Eq. (2). To obtain, for example. A, let k Q = k^ = 

C z > k 2 = V an 6 d k 3 = S in E q- ( 5 )- 
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As a measure of the productivity of the implicit algorithm to be described 
below compared to that of conventional explicit procedures, the value of the 
Courant number can be computed. It requires the eigenvalues of the matrices 
A, B, and C. They are as follows: 


“1,2,3 * k 0 + uk l + 


VK 2 ♦ f l<3 


( 6 ) 


and 


0 4,5 " k 0 + uk l + vk 2 + r k 3 - C 


+ k * + (t) 


2 1/2 


where for A, k Q = ? t , etc.; for B, k Q = n^, etc.; for C, k Q = 5 ^., etc.; 
and c equals the local speed of sound (c^ = yp/p). 


The Courant number is defined as follows: 



(7) 


where it is assumed that as = an = Ac = 1 and a is the maximum 

max 

of the eigenvalues of all the nodal points. 


For the present problem, Eq. (2(a)) is solved in a time-asymptotic fashion 
with interest only in the steady-state solution. This is a result of the mesh 
rotating with the blade and the fact that a cylindrical coordinate system is 
used. The Q T term of Eq. (2(a)) approaches zero as t becomes large, 
thus establishing a convergence criterion. 

NUMERICAL METHOD 


The numerical algorithm used to solve the conservation-law form of the 
Euler equations is based on a class of completely implicit, noniterative, ADI 
(alternating direction implicit) schemes developed by Lindemuth and Killeen 
(ref. 1), Briley and McDonald (refs. 2 and 3), and Beam and Warming 
(refs. 4 and 5). The particular method is a generalization of a conservative, 
approximate factorization scheme in the "delta" form (ref. 4). The procedure 
has been successfully applied by Steger (ref. 8), Pulliam and Steger (ref. 9), 
and numerous others (refs. 10 and 11) for viscous flows. Use of the implicit 
procedure helps remove the stiffness of the problem introduced by a fine 
mesh. Thus, for this problem the implicit procedure permits an integration 
stepsize large enough to obtain steady-state solutions in considerably fewer 
iterations than conventional explicit procedures. 

As applied to Eq. (2) the implicit, spatially factored algorithm using 
Euler implicit time differencing takes the form 
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( 8 ) 


where A, B, and C are the Jacobian matrices aE/aU, a F /a U » and aG/aU, respec- 
tively, and I is the identity matrix. 65 and ^ are second-order central- 
difference operators, and a and v represent the conventional forward and 
backward difference operators. The quantities multiplied by e -j on the 
left-hand side of the equation represent implicit second-order smoothing terms 
(refs. 8 and 9) while the quantity multipled. by e e on the right-hand side 
represents an explicit fourth-order smoothing term. 

The solution of Eq. ( 8 ) consists of first forming the right-hand term (also 
called the steady-state or explicit part) at each grid point. Each of the im- 

n —1 n 

plicit operators, (for instance II + At <$^A n - e^(J v^J) n }) represent block- 

tridiagonal matrices which must be inverted sequentially to obtain AQ n = qo + 1 _ Qn # 
A block lower-upper-decomposition algorithm, see Isaacson and Keller (ref. 12), 
is used for the inversion process. 

For steady-state calculations AQ n = Q n+ 1 - Q n + 0 and the solution then sat- 
isfies the steady-state finite-difference equations (right-hand side of Eq. ( 8 )). 
Linear stability analysis shows unconditional stability for the "delta" form of 
the implicit approximate factorization algorithm. In actual practice, though, 
time step limitations are encountered, although they are usually much less 
stringent than explicit stability bounds. The smoothing terms have been added 
to control nonlinear instabilities. Choices of the time step At and the 
smoothing coefficients (e e , e^) are usually dictated by experience. Linear ana- 
lysis for the smoothing terms does show that for e -j = 0, e e is bounded (e e < 1/4) 
for stability, but for e -j >>e e , £ e is not constrained to any limit for stability. 

The metric terms are obtained using second-order finite-differences for 
terms such as, z^, in Eqs. (4). Fourth-order accuracy in the steady-state 
can be obtained in an efficient manner by introducing fourth-order finite dif- 
ferences for the convective terms on the right-hand side of Eq. ( 8 ), while 
retaining second-order differences on the left-hand side. For Euler implicit 
time differencing this is a stable and accurate process. For more details on 
the above, see Pulliam and Steger (ref. 9). 
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The major portion of the computational work in an implicit finite differ- 
ence algorithm is contained in the solution of the set of simultaneous equa- 
tions, (i.e., the block tridiagonal equations). When even an implicit algo- 
rithm is applied to a system of partial differential equations, one obtains 
block matrix-vector equations that are complicated and time-consuming to 
solve. A method for uncoupling the solution process through a diagonalization 
of the block-matrix structure has been presented by Pulliam and Chaussee 
(ref. 11). The method which is employed here was originally applied to an 
implicit approximate-factorization algorithm for the twoand threedimensional 
inviscid Euler equations in general curvilinear coordinates. Details on the 
theory can be found in reference 13. 

Pulliam and Chaussee (ref. 11) show that the diagonal algorithm retains 
the stability and convergence properties of the standard algorithm. Also, 
since the right-hand side is identical for the two algorithms, the steady- 
state solution is not affected by the diagonalization. Reductions in CPU time 
of up to 30 percent can be realized by using the diagonal form. 

GEOMETRY AND MESH GENERATION DESCRIPTIONS 

The generation of a computational mesh for calculating the flow through 
and around a propeller configuration consists of three steps. Step one is 
composed of generating an x-y-z system of points given the ccoordinates of the 
nacelle and blades. Step two consists of using the x-y-z locations of step 
one and to form a system of parametric cubic patches. Step three takes the 
patches of step two and, based on the specified clustering and/or stretching 
of points, creates a mesh which is used by the computational code. 

The geometry for the blade, spinner, and nacelle must be specified. The 
geometry is described by the blade twist in degrees, as, relative to the 
blade angle at 75 percent blade radius, the local blade chord ratio, b/D, and 
the leading edge alignment ratio, LEA/D, where D is the tip diameter. The 
thickness ratio distribution for each spanwise section is also required. 

The ultimate goal of the mesh generation scheme is to smoothly transition 
from the blade with its twist and taper to a set of surfaces which have : 1) 

zero twist or taper above an r/D =1, 2) a slope d$/dz = 0 at the stagna- 
tion point of the nacelle, and 3) zero slope at some specified distance down- 
stream from the blade trailing edge. The two transition curves which are 
employed to accomplish this are parabolic and cubic in nature. These curves 
allow for end points and/or slopes at the end points to be specified in the 
transitioning process. 

Once these surfaces have been specified, it is an easy matter to extend 
the surfaces at constant values of any specified outer boundary. Usually this 
boundary is chosen to be far enough removed from the configuration such that 
there are no outer interactions affecting the computation. 

This whole procedure is automated with only outer boundary points and/or 
slopes specified in addition to the coordinates of the blade and nacelle sur- 
faces. 
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Step 2 is a self sufficient part of the grid generation program. An x-y-z 
set of points from step one is all that is necessary to create a series of 
parametric cubic patches which describes the entire mesh system in terms of 
cubic equations. 

In the third and final step, the computational mesh used in the flow field 
code is created from the series of patches from step 2. 

Details of the mesh generation procedures are given in reference 13. 


TANGENCY CONDITIONS 

It is necessary to satisfy the tangency condition at both the blade and 
nacelle surfaces. In the case of the blade surface, this requires the contra- 
variant velocity component W to be zero (no flow through the surface). To 
obtain the cylindrical velocity components velocity components along the 
(x,y,z,t) = constant surface, the following system of equations must be solved 


u 

-1 

(n rWr ),r -l l rVVr l,r 

(5 rW 5 * )/r 


' u V 

V 

= 0 

-<V*-Vz> /r (? zVVz )/r 

-< E zVVz> /r 


V-n t 

w 


(n z ?r _n r c z^ -U z C r -S r e z ) 

U z n r -C r n z ) 


_ W -<t_ 


To solve these equations, a complete set of contravariant velocities, U, 

V, and W, must be known. The velocities u, v and w at the surface, are deter- 
mined by linear extrapolation from known data in the interier. The values of 
U and V are then determined from equations (2(c)) and values of u, v and w are 
calculated from Eqs. (9) with W set to zero. 

Flow velocities at the blade tip and trailing edge are determined in the 
same manner as the rest of the blade. No Kutta condition is applied at these 
locations since the velocity of the fluid relative to the blades can be super- 
sonic. Imposition of periodicity at the mesh point just off the blade effec- 
tively applies the Kutta condition at that location. 

Due to insufficient mesh points to resolve the rapid change from low 
velocity to very high velocity in the leading edge region of the blade, flow 
properties at the leading edge of the blade are determined by averaging the 
values immediately upstream and downstream of the leading edge. 

At the nacelle surface, the contravariant velocity V is required to be 
zero to satisfy the tangency condition. In the same manner as for the blade 
surfaces, U and W are determined on the surface and u, v and w are determined 
from Eq. (2(c)) with V set to zero. At the blade nacelle juncture, both V and 
W are set to zero before determining the surface velocities. 

Values of pressure and density on the blade and nacelle surfaces (except 
for the blade leading edge) are determined by linear extrapolation from values 
in the interior of the computational region. 
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As a check on the accuracy of this procedure the fluid total rothalpy is 

monitored throughout the flow field. The total rothalpy H is defined as the 

r 

total enthalpy in the rotating frame of reference or: 


H = — r £ + i (u 2 + v 2 + w 2 - 2r u w) 

r y - i p c 


(10) 


In a flow with uniform free stream, the total rothalpy is constant throughout 
the flow field. 


SUBSONIC OUTFLOW 


The relationships enforcing boundary conditions at the subsonic outflow 
are determined from the governing partial differential equations using the 
method of characteristics. The equations, in non-conservation-law form cast 
in cylindrical coordinates for nonviscous, non-heat conducting, adiabatic, 
compressible flow, under the transformation of Eq. (1) are: 

Continuity 

iE + u +. V IP + w IE + pC 2 \b H + H 3£ + 13W31 + 3U3n 
3t at an . at L 3 ? Mz at Mr r at a<i> an az 
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8 V 

an + 1 

aw 

an + au 

at + av 

at + 1 

aw 

at 


ar r 

an 

a«t> at 

az at 

ar r 


a<t> 


= 0 


( 11 ) 


t- Momentum 

iii+lliE+v — + u — +1 |lPiL+l£ill+lP ii.] = o 
at u at an at p |_a t az an az at azj 


( 12 ) 


n- Momentum 


l^+ll 1Z+V- + W — +- TiP H + l£lH + l£lll = o 
ax at an a? p at ar an ar at ar 


(13) 


t-Momentum 


aw + II w + V in + u — + - - Tie il + IE il + ap i£ 
3t u at v an at p r at a<t> an a<t> at 34* 


(14) 


where the contravariant velocities are defined in Eq. (2(c)). The details for 
deriving the characteristics and the compatibility equations from Eqs. (11 to 14) 
can be found in reference 13. 
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For the outflow (exit) boundary conditions, which is a finite distance 
downstream of the blade trailing edge, a constant pressure is specified. The 
U-c characteristic brings information into the control volume from outside 
when U < c. Thus one dependent variable such as the back pressure must be 
specified. The remaining dependent variables are determined by using the up 
running characteristic (+) in the compatibility relationship to solve for the 
axial velocity, u. The momentum equations are used to solve for the remaining 
velocity components, v and w, and the energy equation is solved for the 
density. 

The uprunning characteristic form of the compatibility equation, with the 
appropriate x and v^'s is written as 



z r \ r 
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where 


R 1 = " v Pn 
R 2 “ - ?U n 
R 3 = - Vv n 



With p specified as a constant in time, p T = 0, Eq. (15) is solved for u T 



For the grid being used for this problem, Cp = = 0 at the down-stream boundary, 

causing a simplification of Eq. (16). If tnis was not the case, Eq. (16) would 
have to be solved simultaneously with the v-and w-momentum equation. 

Writing equations with the simplification, 5 r = £<t> = 0, the following 
equation for the u-velocity is obtained. 


R 1 5 z Pc 

u T - - ( u + C SK + R 2 + TJ p 


(17) 


13 




INITIAL CONDITIONS 


The solution procedure can be started from specified freestream 
conditions. At time t = 0, the values of the dependent variables are set 
equal to freestream conditions at every node of the mesh. At the first time 
step tangency conditions are instantaneously applied in the manner previously 
described. No startup problems have been encountered with this approach. 

An alternate approach is to start from a previously calculated solution. 

If different free stream conditions are specified, the upstream and radial far 
field conditions are changed and the solution procedure is restarted. If a 
different mesh is being used, the dependent variable vector Q must be rescaled 
by the new Jacobian to avoid start up problems. 


NUMERICAL RESULTS 

This analysis has been used to calculate the flow around propellers oper- 
ating at high subsonic speeds. Results will be shown for the eight bladed 
SR-3 propeller (ref. 14) shown in figure 2. The grid used for these calcula- 
tions is shown in figure 3. Figure 3(a) shows the grid along one surface of a 
blade (shaded region) and a portion of the grid in the vicinity of the blade 
at a constant value of the index in the circumferential direction. 

Figure 3(b) shows the grid between two blades at a constant value of the index 
in the radial direction and in figure 3(c) at a constant value of the index in 
the axial direction. A total of 45 points in the "axial" direction, 21 points 
in the "radial" direction and 11 points between adjacent blades in the 
"circumferential" direction were used for the calculations to be described. 

On the blade 21 points in the chordwise direction and 12 points in the 
spanwise direction were used. 

Blade geometry used in the calculations includes deformation due to cen- 
trifugal loads during operation. Blade deformation was calculated using a 
finite element structural analysis program at the design condition of the 
propeller: free stream Mach number of 0.8, advance ratio of 3.06 and blade 

angle at three-quarter radius of 57°. The same deflection was used for all 
blade angles at the design rotational speed. At other rotational speeds the 
deflection was scaled by the square of the rotational speed. 

Flow properties on the blade surfaces of the SR-3 propeller at a free 
stream Mach number of 0.6 are shown in figures 4 and 5. The blade angle at a 
radial location equal to three-quarters of the tip radius is 60.5° and 
the advance ratio is 3.06. Isobars on the suction and pressure surfaces of 
the blades are shown in figure 4. Regions of constant relative Mach number 
are shown in figure 5. The generally increasing level of relative Mach number 
from the nacelle surface toward the tip is a result of the increasing rota- 
tional speed with increasing radius. Each shade of gray indicates a range of 
0.1 in both static pressure ratio and relative Mach number. Both the static 
pressure, figure 4(a), and relative Mach number distribution, figure 5(a), on 
the suction side of the blade indicate higher load on the forward portion of 
the blade. Near the nacelle surface on the suction side, the high pressure 
and low Mach number regions at about the quarter chord location are caused by 
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the rapid increase in spinner diameter at that location. On the pressure side 
of the blade, figure 5(b), the high relative Mach numbers at about mid-chord, 
which is the location of maximum thickness on the blades, is apparently 
related to the thickness distribution of the airfoils used in the blade. 

Computed values of swirl angle downstream of the propeller are shown in 
figure 6. Also shown are experimentally measured values of swirl angle ob- 
tained during the wind tunnel tests described in reference 14. These values 
were measured with an instrumented wedge mounted on a translating probe 
figure 2. These probe results have not been previously reported. Both the 
computed and measured values correspond to an axial location 0.21 diameters 
downstream of the pitch change axis. Although the level of the predicted 
results is considerably higher than the experimental results, the radial 
variation of swirl is in reasonable agreement between the two sets of results. 
An estimate of viscous effects has been made using the compressor methodology 
of reference 15. At a radius ratio r/R = 0.52 which corresponds to a solidity 
(chord/gap) of 1.0, the deviation angle was determined to be 3.15°. This 
corresponds to a reduction of the swirl angle to a value of 5.8° from the 
originally predicted 8.9°. It thus appears that at lower free stream Mach 
numbers the overprediction of power coefficient is due to neglecting viscous 
effects in the analysis and that boundary layer growth on the blade surface 
causes a reduction in swirl angle downstream of the propeller. The unusual 
behavior of the data near r/R=1.0 is apparently caused by the tip vortex roll 
up. The analysis has not predicted this feature of the flow because the mesh 
is too coarse in this region. 

A comparison of computed and measured propeller power coefficient is shown 
in figure 7. Experimental results are taken from reference 14. The power 
coefficient is considerably overpredicted. Since power coefficient is closely 
related to swirl angle the overprediction of power coefficient is consistent 
with the overprediction of swirl angle shown in figure 6. 

Predicted blade surface flow properties for the SR-3 propeller at a free 
stream Mach number of 0.8 are shown in figures 7 and 8. The blade angle at 
three-quarter radius is 61.3° and the advance ratio is 3.06. The pressure 
distribution on the suction side of the blade is shown in figure 8(a). Near 
the nacelle surface the low pressure region at about one quarter chord is 
caused by the rapid increase in spinner diameter at that location. Just down- 
stream a second lower pressure region is observed followed by a rapid change 
to high pressure suggests the presence of a shock wave at about two-thirds 
chord. With increasing distance from nacelle to tip the strength of this com- 
pression decreases and moves to the trailing edge of the blade. On the pres- 
sure side of the blade, figure 8(b), the rapid compression near the nacelle 
surface at about two-thirds chord moves forward with increasing distance from 
nacelle to tip. 

The relative Mach number distribution on the suction side of the blade, 
figure 9(a), generally follows the same trends as the pressure distribution of 
figure 8(a). The level of Mach number near the nacelle surface at about 
two-thirds chord indicates the presence of a shock wave at that location. On 
the pressure side of the blade, figure 9(b), a shock wave is indicated at 
about the same location. Whereas the lines of constant pressure and constant 
relative Mach number on the suction side of the blade have similar shapes, 
they are somewhat different on the pressure side. The reason for this dif- 
ference is unknown at this time. 
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The static pressure distribution in the flow field between two adjacent 
blades is shown in figure 10. The suction surface is at the bottom of the 
plots and the pressure surface at the top. At the nacelle surface, figure 
10(a), the shock wave seen in figures 8 and 9 on both sides of the blades can 
clearly be seen spanning the passage between the blades. A small expansion 
region can be seen near the leading edge of the blade on the suction surface. 
The static pressure distribution on a mesh surface which has r/R=0.82 at the 
mid chord is shown in figure 10(b). The shock wave from the trailing edge of 
the blade on the suction surface spans the blade passage and intersects the 
pressure side of the adjacent blade. 

Relative Mach numbers in the flow field between two blades at a free 
stream Mach number of 0.8 are shown in figure 11. General features are simi- 
lar to the pressure contours shown in figure 10. The main difference is near 
the pressure surface in figure 11(b) where a rise in Mach number is observed 
without a corresponding drop in pressure in figure 10(b). The reason for this 
is not known. The results shown in figure 11(b) are in very good qualitative 
agreement with laser velocimeter measurements presented in reference 16. A 
significant difference is the location of the shock wave. In the computed 
results the shock wave originates at the trailing edge of the blade on the 
suction surface whereas the data of reference 16 indicates a shock wave loca- 
tion somewhat upstream of the trailing edge. This behavior is consistent with 
shock wave boundary layer interactions in other types of flow fields in which 
viscous effects cause the actual shock wave location to be upstream of the 
location predicted by inviscid analyses. 

Predicted and measured swirl angle downstream of the propeller at a free 
stream Mach number of 0.8 are shown in figure 12. Numerical results for swirl 
angle are higher than those measured in the wind tunnel (ref. 14) similar to 
the results at a Mach number of 0.6 (fig. 6). The observed shock wave loca- 
tion upstream of the trailing edge would result in reduced blade loading 
causing reduced swirl relative to a trailing edge shock wave location. This 
would at least partially explain the overprediction of swirl angle shown in 
figure 12. This overprediction is approximately equal to five degrees whereas 
at a free stream Mach number of 0.6 the discrepancy is about three degrees. 
This implies that at free stream Mach number of 0.8 a mechanism in addition to 
boundary layer growth is causing decreased flow turning. This mechanism could 
be the shock wave boundary layer interaction discussed previously. The com- 
pressor methodology of reference 15 cannot be used to estimate viscous losses 
at this condition. This methodology is based on low speed cascade data and 
the presence of shock waves in the flow field clearly makes it not applicable. 
The radial variation of swirl angle agrees reasonably well with the experi- 
mental results. The region with the largest discrepancy is again the tip 
region where the coarse mesh cannot resolve the details of the flow associated 
with the tip vortex. 

The general features of the tip vortex, however, are predicted by the 
analysis as shown in figure 13. The secondary flow velocity vectors on a mesh 
surface just downstream of the propeller at a free stream Mach number of 0.8 
as predicted by the analysis are shown. These vectors are composed of the 
radial and circumferential components of the flow field velocities. Two blade 
passages are shown with a blade at each edge of the region shown and a third 
blade approximately in the middle. The tip vortex can clearly be seen as well 
as the jump in spanwise component of velocity in the inboard region associated 
with increasing loading with radial distance. 
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A summary of propeller power coefficient for the SR-3 propeller at a free 
stream Mach number of 0.8 is shown in figure 14. The SR-3 propeller was 
designed for a cruise condition at a free stream Mach number of 0.8, advance 
ratio of 3.06 and a power coefficient of 1.7. Power coefficient is shown for 
three blade angle settings. Although the power coefficient is again over- 
predicted, the variation of power coefficient with blade angle is predicted 
quite well. The overprediction of power coefficient is consistent with the 
overprediction of swirl angle shown in figure 12 and would improve if shock 
wave boundary layer interaction effects were included. 

CONCLUSIONS 

A computer program, NASPROP-E, has been written to calculate the flow 
field around high speed propellers. This program has been applied to a pro- 
peller with eight highly swept blades and an axisymmetric nacelle operating at 
transonic conditions. Numerical results indicate the presence of shock waves 
in the flow field and the roll up of the tip vortex. 

Comparisons of numerical results with data indicate the following: 

1. Propeller power coefficient and swirl angle downstream of the pro- 
peller are overpredicted. At lower free stream Mach numbers this discrepancy 
can be attributed to viscous effects. At higher free stream Mach numbers this 
discrepancy is apparently due to shock wave boundary layer interaction effects. 

2. Radial variation of swirl angle downstream of the propeller is pre- 
dicted reasonably well. Details associated with the tip vortex are not as 
well defined due to a relatively coarse grid in that region. 

3. The variation of propeller power coefficient with blade angle is pre- 
dicted very well. 
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Figure 1, ~ Coordinate system. 


Figure 2. - SR-3 propeller in the NASA-Lewis Research Center 8x6 
supersonic wind tunnel. 






(c) f; = constant surface. 
Figure 3. - Concluded. 



(a) Suction side 


(b) Pressure side 


Figure 4. - Predicted static pressure distribution, p/p,*,, on both sides of the SR-3 propeller blade. Free stream Mach number, 0.6; advance ratio, 3.06; 
blade angle at three-quarter radius, 60.5°. 




(a) Suction side. 


(b) Pressure side. 


Figure 5. - Predicted distribution of relative Mach number on both sides of the SR-3 propeller blade. Free stream Mach number, 0.6; advance ratio, 3.06; 
blade angle at three-quarter radius, 60.5° . 





Figure 6. - Comparison of predicted and mea- 
sured swirl angle downstream of the SR-3 
propeller blade. Axial location, 0.21 diam- 
eters downstream of the pitch change axis; 
free stream Mach number, 0.6; advance 
ratio, 3.06; blade angle at three-quarter 
radius, 60.5°. 
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Figure 7. - Comparison of predicted and measured 
power coefficient for the SR-3 propeller. Free stream, 
Mach number, 0.6; advance ratio, 3.06; blade angle 
at three-quarter radius, 60.5°. 








(a) Suction side 


(b) Pressure side 


Figure 9. - Predicted distribution of relative Mach number on both sides of the SR-3 propeller blade. Free stream Mach number, 0.8; advance ratio, 3.06; 
blade angle at three-quarter radius, 61.3° . 
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Figure 10. - Predicted static pressure distribution, p/p^, between the blades of the SR-3 propeller. Free stream Mach number, 0.8; advance ratio, 3.06; 
blade angle at three-quarter radius, 61.3°. 
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(a) r/R = 0.24 (nacelle surface). 


(b) r/R = 0.82. 


Figure 11. - Predicted distribution of relative Mach number between the blades of the SR-3 propeller. Free stream Mach number, 0.8; advance ratio, 
3.06; blade angle at three-quarter radius, 61.3° . 
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Figure 12. - Comparison of predicted and mea- 
sured swirl angle downstream of the SR-3 
propeller blade. Axial location, 0.21 diam- 
eters downstream of the pitch change axis-, 
free stream Mach number, 0.8; advance 
ratio, 3.06; blade angle at three-quarter 
radius, 60.5°. 



Figure 13. - Predicted secondary flow velocity vectors down- 
stream of the SR-3 propeller. Free stream Mach number, 
0.80; advance ratio, 3.06; blade angle at three-quarter 
radius, 61. 
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Figure 14 - Comparison of predicted and measured 
power coefficient for the SR-3 propeller. Free stream 
Mach number, 0.8 ; advance ratio, 3.06. 
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